/*
Network Design Project
"Minimum Bound Degree Spanning Tree"
Version: 0.4
Authors: Alessandro Baroni, Paolo Bozzelli, Giuseppe Cascavilla
Domain: univaq.it
*/

#define USAGE "mbdst-formulation.out <min_n> <max_n> <min_w> <max_w> <k>\n"

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <time.h>
#include <glpk.h>
#include "lib/ourlib.h"

int next_comb(int[], int, int);
void process_c(int[], int);

glp_prob *lp;
glp_graph *G;
int e;
int n;

int main(int argc,char *argv[]) {
	
	const char* LP_FILE = "mbdst-formulation.lp";
	const char* GRAPH_FILE = "generated_graph.txt";
	int i,k,final_solution_arcs;
	int *comb;
	double *x,z;
	
	lp = glp_create_prob();
	G = custom_create_graph();
	i = glp_read_lp(lp,NULL,LP_FILE);
	i = custom_read_graph(G, GRAPH_FILE);
	e = G->na;
	n = G->nv;
	
	
	glp_create_index(lp);
	
	/* 
	*	k is the size of the subsets; for {1, 2}, {1, 3}, ... it's 2 
	*	n is the size of the set; for {1, 2, 3, 4} it's 4
	*/
	

	comb = (int*)malloc(sizeof(int)*(n+1));	/* comb[i] is the index of the i-th element in the combination */
	for(k = 2;k<=n;k++){
		int i;
		/* Setup comb for the initial combination */
		for (i = 0; i < k; ++i){
			comb[i] = i+1;
		}
		
		/* Print the first combination */
		process_c(comb, k);

		/* Generate and print all the other combinations */
		while (next_comb(comb, k, n)){
			process_c(comb, k);
		}
	}
	//TODO LIBERARE LA MEMORIA!!!
	//free(comb);
	
	glp_write_lp(lp, NULL, "mbdst-formulation-ex.lp");
	glp_simplex(lp, NULL);
	
	final_solution_arcs = 0;
	x = malloc (sizeof(double)*glp_get_num_cols(lp));
	z = glp_get_obj_val(lp);
	for (i = 0; i < glp_get_num_cols(lp); i++) {
		x[i] = glp_get_col_prim(lp, i+1);
		if (x[i] != 0.0) {
			final_solution_arcs++;
		}
	}
	
	// Stampa su stdout del risultato
	printf("\nObjective:\tz = %g;\n", z);
	printf("Vector:\n");
	for (i = 0; i < glp_get_num_cols(lp); i++) {
		//printf("\tx%d = %g;\n", i+1, x[i]);
		printf("\t\t%s = %g\n", glp_get_col_name(lp, i+1), x[i]);
	}
	
	printf("La soluzione finale ha %d archi.\n", final_solution_arcs);
	custom_write_sol(G,lp,"soluzione");
	
	free(x);
	/*
	printf("free\n");
	glp_delete_prob(lp);
	printf("free\n");
	custom_delete_graph(G);
	printf("free\n");
	*/
	return 0;

}


/* Prints out a combination like {1, 2} */
void process_c(int comb[], int k) {
	int row_num,i,j;
	int *this_row_indexes;
	double *this_row_coeffs;
	int this_row_len;
	glp_arc* arcp;
	
	this_row_len = 0;
	this_row_indexes = malloc(sizeof(int)*(e+1));
	this_row_coeffs = malloc(sizeof(double)*(e+1));
	
	/*
	for(i = 0; i<k;i++){
		printf("i=%d\t comb[i]=%d\n",i,comb[i]);
	}
	*/		
	for(i = 0; i<k;i++){
		// per ogni nodo in comb
		arcp = G->v[comb[i]]->out;
		//printf("processing node %d\n",comb[i]);
		while (arcp != NULL) {	
			int v1,v2,v1_is_valid,v2_is_valid;
			//per ogni arco uscente dal nodo corrente
			v1 = arcp->tail->i;
			v2 = arcp->head->i;
			v1_is_valid = 0;
			v2_is_valid = 0;
			//printf("processing (%d,%d)\n",v1,v2);
			//se (i,j) ha i e j in comb
			for(j = 0; j<k;j++){
				if(v1 == comb[j]) v1_is_valid = 1;
				if(v2 == comb[j]) v2_is_valid = 1;
			}
			
			if(v1_is_valid && v2_is_valid){
				//genero la stringa "x_i_j"
				char arcp_var[128];
				char tmp[128];
				int var_index;
				strcpy(arcp_var,"x_");
				itoa(v1,tmp);
				strcat(arcp_var,tmp);
				strcat(arcp_var,"_");
				itoa(v2,tmp);
				strcat(arcp_var,tmp);
				//ripesco il numero della variabile "x_i_j" e lo metto in var_index
				var_index = glp_find_col(lp,arcp_var);
				this_row_len++;
				this_row_indexes[this_row_len] = var_index;
				this_row_coeffs[this_row_len] = 1.0;
				//printf("%s is %d\n",arcp_var,var_index);
			}
			arcp = arcp->t_next;
			//end while
		}
		//end for
	}
	
	//aggiungo la riga ad lp
	glp_add_rows(lp,1);
	row_num = glp_get_num_rows(lp);
	glp_set_mat_row(lp, row_num,this_row_len,this_row_indexes,this_row_coeffs);
	//glp_set_row_name(lp, rownum, constraint_name);
	glp_set_row_bnds(lp, row_num, GLP_UP, 0.0, k-1);	// <= |S|-1
	
	free(this_row_indexes);
	return;
}

/*
	next_comb(int comb[], int k, int n)
		Generates the next combination of n elements as k after comb

	comb => the previous combination ( use (0, 1, 2, ..., k) for first)
	k => the size of the subsets to generate
	n => the size of the original set

	Returns: 1 if a valid combination was found
		0, otherwise
*/
int next_comb(int comb[], int k, int n) {
	int i = k - 1;
	++comb[i];
	while ((i >= 0) && (comb[i] >= n - k + 1 + i)) {
		--i;
		++comb[i];
	}

	if (comb[0] > n - k) /* Combination (n-k, n-k+1, ..., n) reached */
		return 0; /* No more combinations can be generated */

	/* comb now looks like (..., x, n, n, n, ..., n).
	Turn it into (..., x, x + 1, x + 2, ...) */
	for (i = i + 1; i < k; ++i)
		comb[i] = comb[i - 1] + 1;

	return 1;
}
